library(nycflights13)
library(mosaicData)
library(janitor)
library(tidyr)
library(broom)
library(mdsr)
library(RColorBrewer)
library(osmdata)
library(ggmap)
library(zipcodeR)
library(tidyverse)
library(paletteer)
library(gt)

# Find data with na's
cancelled_flights <-
    flights %>% 
    filter(is.na(flights$dep_time) & is.na(flights$arr_time) & 
               is.na(flights$air_time)) %>% 
    select(month,day,arr_time,air_time,dep_time)

#Find flights cancelled
number_flights_cancelled <- 
  cancelled_flights %>% 
  group_by(month) %>% 
  summarise(n_flights_cancelled = n())

# Find total_number_flights
total_number_flights <- 
  flights %>% 
  group_by(year,month) %>% 
  summarise(n_flights_total = n())

# Join tables
proportion_flights_cancelled <-
  inner_join(number_flights_cancelled,total_number_flights,
             by = c("month"))

proportion_flights_cancelled %>%
    mutate(proportion_flight_cancelled = 
                 ((n_flights_cancelled/n_flights_total))) %>% 
    round(2) %>% 
    arrange(desc(proportion_flight_cancelled)) %>% 
    rename(Month = "month",flights_cancelled = "n_flights_cancelled",
           Total_flights = "n_flights_total",
           Proportion_flights_cancelled = "proportion_flight_cancelled") %>%
    select(year,everything()) %>% 
    gt() %>% 
    fmt_percent(Proportion_flights_cancelled, decimals = 0)

weather %>% 
  group_by(year,month,day) %>% 
  mutate(
        origin = as.factor(origin),
        year = as.factor(year),
        month = as.factor(month),
        day = as.factor(day),
        hour = as.factor(hour)) %>% 
    summarize(across(where(is.numeric), mean , na.rm = T)) %>% 
    filter(precip > 0) %>% 
    group_by(month) %>% 
    summarise(days_with_rain = n()) %>% 
    arrange(desc(days_with_rain)) %>% 
    rename(Month = "month",Days_with_rain = "days_with_rain") %>% 
    gt() 

flights %>%
  select(-year) %>% 
      inner_join(planes, by = c('tailnum' = 'tailnum')) %>% 
      select(tailnum,year, everything()) %>% 
      group_by(tailnum,year) %>% 
      summarise(Mean_airtime = mean(air_time)) %>%
      arrange(year) %>% 
      drop_na() %>% 
      head(5)

# Load data
data("HELPfull")

HELPfull %>% 
    clean_names()  %>% 
    filter(id == 1) %>% 
    filter(time %in% c(0,6)) %>% 
    pivot_longer(cols = c(drugrisk,sexrisk), 
                 names_to = "type_of_risk", 
                 values_to = "risk_score") %>% 
    select(id,type_of_risk,risk_score,time,everything()) %>% 
    arrange(type_of_risk) %>% 
    select(1:9) %>% 
    gt()

table <- 
HELPfull %>%
    clean_names() %>% 
    select(ces_d,age,secd_sub) %>% 
    na.omit() %>% 
    mutate(secd_sub = as.factor(secd_sub)) %>%
    nest(data = -secd_sub) %>% 
    mutate(
      fit = map(data, ~ lm(ces_d ~ age, data = .x)),
      tidied = map(fit, tidy, conf.int =T)) %>% 
    unnest(tidied) %>% 
    select(-c(data,fit)) %>% 
    filter(term %in% "age") %>% 
    clean_names() %>% 
    arrange(secd_sub) %>% 
    drop_na() 

levels(table$secd_sub) <- c("None","Alcohol","Cocaine" ,
                            "Heroine","Barbituates", "Benzos", 
                             "Marijuana", "Methadone", "Opiates")

table %>% gt()

link <- "https://raw.githubusercontent.com/opencasestudies/ocs-healthexpenditure/master/data/KFF/healthcare-coverage.csv"

data_csv <-  
    read.csv(file = link, skip = 2, nrows = 52, 
             na.strings = "N/A")  %>% 
    clean_names() 
    

# Remove x from colnames
colnames(data_csv) <- sub("^x", "", colnames(data_csv))

data_csv <- 
  data_csv %>%   
      mutate(across(starts_with("201"), as.integer))

data_csv <- 
    data_csv %>%
        pivot_longer(-location, names_to = "year", 
                     values_to = "tot_coverage") %>% 
        separate(year,c("year","type"),extra = "merge") %>% 
        mutate(year = factor(year),
               type = factor(type),
               type = str_to_sentence(type))

data_csv$type <-  sub("_", "-", data_csv$type)

head(data_csv, 10) %>% 
  gt()

# Load data
data("ViolationCodes")
data("Violations")

# Join ViolationCodes and Violations datasets

data_nyc_restaurants <- 
    inner_join(ViolationCodes,Violations, 
             by = c("violation_code" = "violation_code")) %>% 
    
    # Clean dataset delete columns
    select(-c(violation_description,cuisine_code,phone,action,
              street,grade_date,boro,record_date,
              grade,inspection_type,violation_code,building))  %>% 
    
    drop_na() %>% 
    mutate(across(where(is.character), as.factor),
           camis = factor(camis)) %>% 
    select(zipcode,camis,dba,inspection_date,everything()) %>% 
    arrange(zipcode,camis,inspection_date)  

data_per_rest <-
        data_nyc_restaurants %>% 
        drop_na() %>% 
        group_by(zipcode,camis,inspection_date) %>%
        
        #Create critical flags column
        mutate(n_critical_flags = 
                   ifelse(critical_flag == "Critical",TRUE,FALSE)) %>%
        filter(n_critical_flags > 0) %>% 
    
        summarise(number_of_critical_flags = 
                    sum(n_critical_flags == 'TRUE'),
                score = mean(score)) %>% 
      
        arrange(zipcode,camis)

data_number_visits_per_rest <- 
    data_per_rest %>%
        # How many times a restaurant was visit?
        group_by(zipcode,camis) %>%
        summarise(number_visits = n(),
                  mean_score = mean(score),
                  total_number_of_critical_flags = 
                      sum(number_of_critical_flags)) %>% 
  
        arrange(zipcode,camis) %>% 
        select(zipcode,camis,
               total_number_of_critical_flags,everything())

data_agregated_by_zipcode <-
    data_number_visits_per_rest %>% 
    #select(-camis) %>% 
    group_by(zipcode) %>% 
    summarise(mean_inspector_visits = mean(number_visits),
              mean_inspector_score = mean(mean_score),
              mean_critical_violations = 
                  mean(total_number_of_critical_flags) ) %>% 
    arrange(zipcode)  

data_agregated_by_zipcode %>% 
    head(10)

ggplot(data = data_agregated_by_zipcode ,
       aes(x = mean_inspector_visits, y = mean_critical_violations)) +
    geom_point() +
    geom_smooth(method = lm) +
    theme(axis.text.y  = element_text(size = 14),
                # Legend position and Axis size 

        axis.text.x   = element_text(size = 14),
        axis.title.y  = element_text(size = 14),
        axis.title.x  = element_text(size = 14),
        # Add borders to the plot
        panel.border = element_rect(colour = "black", 
                                    fill= NA,size = 1.3)) +
    ylab("Critical Violations") + xlab("Inspector Visits")

# Get lat long for each zipcode
zip_long_lat <- geocode_zip(data_agregated_by_zipcode$zipcode) %>% 
        mutate(zipcode = as.integer(zipcode))

# Join agregated data and lat long data
nyc <- 
    inner_join(data_agregated_by_zipcode,zip_long_lat, 
                  by = c("zipcode" = "zipcode")) %>%
    select(zipcode, lat,lng, everything()) %>% 
    drop_na()

# Get map and set palette
nyc_map <- get_map(getbb("New York"), maptype = "toner-background")
my_palette <- colorRampPalette(rev(brewer.pal(11, "Spectral")))

ggmap(nyc_map) +
  geom_point(data = nyc, aes(x = lng, y = lat, 
                             color = mean_critical_violations),
             size = 9) +
  scale_colour_gradientn(colours = my_palette(10), limits=c(1, 22)) +
  theme(legend.position = "bottom",
        legend.key.size = unit(2, 'cm'),
        legend.title = element_text(size=30),
        legend.text = element_text(size=30)) +
    labs(color = 'Critical \n Violations')